↜ Back to index Introduction to Numerical Analysis 1
Solution Lecture a4
Example solutions of Report 1.
Exercise 1
! Implementation of the midpoint method for y' = y sin t
implicit none
integer i, n
real y, h, t, f, k
n = 20
! interval (0, 5)
h = 5. / n
! Initial data
y = 1.
print *, 0., y
do i = 1, n
! One step of the midpoint method.
t = (i - 1) * h
f = y * sin(t)
k = y + 0.5 * h * f
t = (i - 0.5) * h
f = k * sin(t)
y = y + h * f
print *, i * h, y
enddo
endExercise 2
! Implementation of the backward Euler method for y' = -10(y - cos t)
implicit none
integer i, n
real y, h, t, f
n = 10
h = 5. / n
! Initial data
y = 0.
print *, 0., y
do i = 1, n
t = i * h
y = (y + 10 * h * cos(t)) / (1 + 10 * h)
print *, t, y
enddo
endExercise 3
1-2 without arrays
! Implementation of the midpoint method for a system y_1' = y_2, y_2' = - y_1
implicit none
integer i, n
real h, t
real y1, y2, k1, k2
n = 20
h = 10. / n
! Initial data
y1 = 1.
y2 = 0.
print *, 0., y1, y2
do i = 1, n
k1 = y1 + 0.5 * h * y2
k2 = y2 + 0.5 * h * (-y1)
y1 = y1 + h * k2
y2 = y2 + h * (-k1)
print *, i * h, y1, y2
enddo
end1-2 with arrays
! Implementation of the midpoint method for a system y_1' = y_2, y_2' = - y_1
implicit none
integer i, n
real h, t
real y(2), k(2)
n = 20
h = 10. / n
! Initial data
y = [1., 0.]
print *, 0., y
do i = 1, n
k = y + 0.5 * h * [y(2), -y(1)]
y = y + h * [k(2), -k(1)]
print *, i * h, y
enddo
end1-3
! Implementation of the midpoint method for the system y_1' = y_2, y_2' = - y_1
implicit none
integer i, n, m
real h, t
real y(2), k(2)
do m = 3, 10
n = 2**m
h = 10. / n
! Initial data
y = [1., 0.]
do i = 1, n
k = y + 0.5 * h * [y(2), -y(1)]
y = y + h * [k(2), -k(1)]
enddo
print *, n, sqrt((y(1) - cos(10.))**2 + (y(2) + sin(10.))**2)
enddo
endExercise 4
! Implementation of the backward Euler method for the system y_1' = y_2, y_2' = - y_1
implicit none
integer i, n
real h, t
real y(2)
n = 100
h = 10. / n
! Initial data
y = [1., 0.]
print *, 0., y
do i = 1, n
y = (y + h * [y(2), -y(1)]) / (1 + h**2)
print *, i * h, y
enddo
end